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^ , Abstract 

^\| , The Fokker-Planck equation is considered, which is connected to the birth and 

death process with immigration by the Poisson transform. The fractional derivative 

\l \ in time variable is introduced into the Fokker-Planck equation. From its solution (the 

J^ ■ probability density function), the generating function (GF) for the corresponding 

^— V \ probability distribution is derived. We consider the case when the GF reduces to that 

^T^ ■ of the negative binomial distribution (NBD), if the fractional derivative is replaced 

0^ \ to the ordinary one. Formulas of the factorial moment and the Hj moment are 

^D • derived from the GF. The Hj moment derived from the GF of the NBD decreases 

^^ \ monotonously as the rank j increases. However, the Hj moment derived in our 

^^ ■ approach oscillates, which is contrasted with the case of the NBD. Calculated Hj 

i-Q \ moments are compared with those given from the data in pp collisions, e^e~ and ep 
collisions. A phenomenological meaning of introducing the fractional derivative in 

^ i' time variable is discussed. 
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1 Introduction 



The negative binomial distribution is often used for the analysis of observed multiplicity 
Cd ' distributions in high energy hadron-hadron {hh) and e^e^ collisions. The cumulant mo- 

ment derived from the generating function of the negative binomial distribution (NBD) 
does not show oscillatory behaviors as the rank of the cumulant moment increases. On 
the other hand, cumulant moments obtained from observed multiplicity distributions in 
hh and e+e^ collisions show oscillatory behaviors |]l|, |2|. Those behaviors can be ex- 
plained if multiplicity distributions truncated at the highest observed multiplicities are 
used for the calculation of Hj moments: In hh collisions, calculated results from the 
NBD and those from the modified NBD both fit the data well [^, ^]. In e^e~ collisions, 
calculated Hj moments by the use of the modified NBD describe the oscillatory behavior 
of the data well. However, those by the NBD oscillate much weaker than the data, and 
cannot explain the behavior of the data |^] . 

The NBD and the MNBD are derived from the branching equations; the former is 
from a birth and death process with immigration, and the latter is from a pure birth (or 



birth and death process). In those branching equations, it is assumed that particles are 
produced instantaneously, in other words, without memory. 

In high energy particle-particle collision processes, it is considered that a proper time 
is needed for a secondary produced particle to behave as an independent particle from the 
parent particle after the collision of the parent with a target particle |g . In high energy 
hadron- nucleus or lepton- nucleus collisions, this effect should be observed as a supression 
of multiplicity compared with the case of instantaneous collision, because the incident 
particle can collide with another target particle in the same nucleus within the proper 
time after the first collision. This proper time is called the formation zone [0], which 
means some memory effect should be existent in the high energy particle production 
processes. 

In the branching equations, particles are assumed to be produced successively. If a 
memory effect is introduced into the branching process, it will be very interesting what 
results come out. 

The birth and death process with immigration is described by the following equation, 

dP{n,t) 



dt 



Xo[P{n-l,t)-P{n,t)] 

+A2 [(n - l)P(n - 1, t) - nP{n, t)] 
+Ai [(n + l)P{n + 1, t) - nP{n, t)] 
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where P{n, t) denotes the probability distribution that n particles are existent at time 
t, Ao denotes an immigration rate, Ai a death rate and A2 a birth rate. If the initial 
condition is taken as 

P{n,t = 0) = 0, 

the solution of equation (|l]) becomes the NBD. 

The probability density function (KNO scaling function) ^{z, t) is connected to the 
probability distribution (multiplicity distribution) P{n,t) by the Poisson transform, 



P{n,t) 



z" exp[—{no)z]ilj{z, t)dz. 



(2) 
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The KNO scaling function ip{z,t) is obtained from the multiplicity distribution P{n,t) 
by the inverse Poisson transform. 
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Applying the inverse Poisson transform 
tion, 

dil){z,t) _ d 
dt dz 

where 



to Eq.(|l|), we obtain the Fokker-Planck equa- 
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If time derivative in Eq.(^) is replaced to the fractional one, we have reached to the 
fractional Fokker-Planck equation in time variable as a model for high energy particle 
production processes, in which a memory effect is taken into account. 

The fractional calculus has been investigated for hundreds of years JS], ^. Recently, the 
fractional Fokker-Planck equation in time variable was derived from the continuous time 



random walk |1C]. It is applied to the analysis of anomalous diffusion phenomena |11]. 
The fractional derivative in space variable is introduced into the Fokker-Planck equation 
to describe the Levy process ||T^ . 

We would like to take the fractional Fokker-Planck equation in time variable cor- 
responding to the branching equation (||) as a model for particle production process, 
and to investigate it's solution, which reduces to the gamma distribution when the frac- 
tional derivative is replaced to the ordinary one. We also examine the effect of fractional 
derivative or introducing the memory effect on the behavior of cumulant moments. 

2 The Fractional Fokker-Planck Equation 

In the following we consider the fractional Fokker-Planck equation (FFPE), 

di;{z,t) 
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with the initial condition, 

ij{z,t = 0) = 5iz-zo). (7) 

In Eq.(^), oD( denotes the Riemann-Liouville fractional derivative |^, ^ defined by 

1 d" '•* 



where n is a positive integer. If a = 1, Eq.@ reduces to Eq.(0). 



According to the method proposed by Barkai and Silbey [11|, we assume that 



^{z,t)= / Rs{t)Gs{z)ds, (9) 
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and that function Gs{z) satisfies the following equations, 

LfpGs{z) = §-^Gs{z), 

Go{z) = 6{z-zo). (10) 

Then function Gs{z) is given as, 

G.(z) = -exp[ ][-^^^^^^) h-ii jTp ), (11) 



where 



fc = ^, A = ^, p=l-e-^'. 
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In the fohowing, we assume that 7 = Ai — A2 > 0. Equation dll]) can be expanded as 

(13) 
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where Lm {z) denotes the Laguerre polynomial. 
Applying the Laplace transform to Eq.(0), we find 



uRs{u) + u 
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Gs{z)ds= l-u^-'^Ro{u) 5{z-zo), 



where Rs{u) is the Laplace transform of Rs{t), 
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Rs{u) = / i?,(t)e-"*dt. 
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Furthermore, we assume that each side of Eq.(|l4|) is equal to zero; 

u'-'^Roiu) = 1, _^i-"^Mf) = uRsiu). 

ou 

The solution of Eq.([lq) is given by 

Rs{u) = n°-^exp[-sn"]. 
Then Rs{t), the inverse Laplace transform of Rs(u), is written as 

X fco+ioo 

Rs{t) = TT- / i?s(n)e"*(in 

ZTTI J CQ—ioo 



+— Q rco+ioo 



2m 



CQ—lOO 



a°-^ exp[a - -a"]da, (co > 0). 



Therefore, the solution for the FFPE (0) is given by 
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where Ea{—t) for f > denotes the Mittag-Leffler function [13|, 
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It is written in the infinite series as 
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In the hmit of zq — > +0, Eq.(^) reduces to 
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If a = 1, Eg. ([2^) coincides with the gamma distribution, the KNO scahng function of 
the NBDQ. 

3 Factorial Moment and Cumulant Moment 

The generating function (GF) for the multiphcity distribution P{n, t) is defined as, 

oo 

n(u) = ^P(n,t)M". (23) 

n=0 

The multiphcity distribution and the jth rank factorial moment are given from Eq.(|2 
respectively as 
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= y n(n-l)---(n- j + l)P(n). (24) 
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From Eqs.(|2|) and (p3D, the GF is written by the use of the Laplace transform of the 
KNO scaling function 'ip{z, t) as 

POO 

U{l-u/{nQ))= ^{z,t)e''''dz. (25) 
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Then the GF corresponding to Eq.(|22|) is given as, 

n(„) = t r(^ + A) (-,(„„)(„ -i)r 



^ We have considered the FFPE for 7 = Ai — A2 > 0, where the birth rate A2 is less than the 
immigration rate Ai . In the case for 7 < 0, the solution of the FFPE (0) in the limit of zq -^ +0, is given 
by 



^("'') = f7A)TfcFE^^'"'(m)^"(-(^ + ^)Wn 
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If Q = 1, the above equation coincides with the gamma distribution. However, we cannot calculate the 
factorial moment from it, because the exponential damping factor in z variable is not contained in the 
equation, contrary to Eq.(p2|). 



The multiplicity distribution and the factorial moment are given from Eqs.(^) and ( pq ) 
respectively as, 
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The jth rank normalized factorial moment is given by 

^_ /, _T{X + j)j:L=oi-'^rjC^E^i-m^t^ 
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The fcth rank cumulant moment is defined by the following equation, 
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From Eqs.(|28|), (|29| ) and (30), we obtain the recurrence equation for the Hj moment, 
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where 

^j = '^j/fj- 
If a = 1, Eq.(|2^) reduces to the generating function for the NBD with mean multi- 
plicity (nb) = A;A(no)(l - e"^*), 
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From Eq.(32), the NBD is given as 

r(A + n) ((nb)/A)" 



P{n,t) 



(32) 



(33) 



r(A)r(n + l)(l + (nb)/A)«+^' 
The normalized factorial moment and the Hj moment for the NBD are given respectively 



as, 



r(A + j) r(A + i)(j-i)! 



(34) 
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As can be seen from Eqs.(|29|) and (p4|), difference between the normalized factorial mo- 
ment derived from the FFPE (0 < a < 1) and that of the NBD (a = 1) is given by 
Mittag-Lefher functions. We can see from Eq.([20|) that 



lim Eai-jf") = 0. 

t— »+oo 



Then, equation (p9) coincide with Eq.(34) in the limit of t ^ -|-oo. 



4 Calculated Results 

At first, calculated results of the Mittag-Leffler function Ea{—t) is shown in Fig.l. It 
is a decreasing function of variable i, and as a increases from to 1, it decreases more 
faster as a function of variable t. 

In the following calculations, observed values of (n) and C^(= (n^)/(n)^) for the 
charged particles are used. Then, if a and 7f" are given, A in Eq.(p^) is determined by 
the following equation. 



In order to see the effect of the fractional derivative, i.e. < a < 1, to the oscillatory 
behavior of Hj moments, calculated Hj moments are shown in Fig. 2a with a =0.25, and 
7t° =1.5, 2.0 and 2.5, in Fig.2b with a =0.50, and 7*" =1.0, 1.5 and 2.0, and in Fig.2c 
with a =0.75, and 7^ =0.5, 1.0 and 1.5. In our calculation, (n) = 29.2 and C^ = 1.274, 
observed values in pp collisions at ^/s = 546 GeV, are used [14|. If q is fixed, oscillation 
of the Hj moment as a function of rank j becomes weaker as the value of parameter 7^" 
increases. If 7t" is fixed, oscillation of Hj moments become much weaker as a increases 
from to 1. 

In Fig. 3, our calculation with a = 0.5 and 7^^ = 2.23 is compared with the Hj 



moment obtained from the data in pp collisions at y/s = 546GeV [14|. Parameter 7^°^ 
is adjusted with a step of 0.01 so that the first relative minimum of the calculated Hj 
moment should be located near the rank of that obtained from the data as much as 
possible. The calculated first relative minimum value is Hj = —3.32 x 10~^, and the 
absolute value of it is much smaller than that obtained from the data. However, we can 
see from Figs. 2b and 3, the calculated Hj moment with a = 0.5 and 7t" = 1.0 oscillates 
as strong as that from the data. 

In Fig. 4, the calculated Hj moment with a = 0.5 and 7t" = 90.0 is compared with 
that obtained from the data in e+e~ collisions at -y/s = 91 GeV |15|. As can be seen 



from the figure, the calculated value of the first relative minimum is almost the same 
with the data, and the strength of the oscillation of calculated Hj moment is comparable 
with the data. 



The Hj moment in ep collisions in the interval 185 < W < 220 GeV |16] is also 
analyzed, and the results are shown in Fig5. Calculated Hj moment with a =0.5 and 
7t" =11.2 (A =20.00) well reproduces the first relative minimum of the data. For com- 
parison, the Hj moment calculated with the NBD ( |33| ) truncated at the highest observed 
charged multiplicity is also shown. The parameters of the truncated NBD are deter- 
mined by the minimum chi-squared fit with the observed charged multiplicity distribu- 
tion. The first relative minimum of the calculated Hj moment with the truncated NBD is 
H^ = —3.42 X 10~^ which is different from the first relative minimum i^s = —4.07 x 10~^ 
obtained from the experimental data. However the strength of the oscillation of Hj mo- 
ments calculated with the truncated NBD is comparable with the data. Estimation of 
Hj moment is obtained by the use of Eq.(pl|) from the Fj moment both in the theoret- 
ical calculation and in the calculation from the experimental data. In order to see the 
relation between the behavior of the Hj moment and that of the normalized factorial 



moment Fj as a function of rank j, we also analyze the Fj moment in ep collisions in 



the interval 185 < T^ < 220 GeV [If]. The parameters are the same with those in Fig. 5. 



The calculated Fj moment with Eq.(15) is compared with the data in Fig. 6a. Our cal- 
culated result in the normalized factorial moment well reproduce the experimental data 
up to the fourth rank. The difference between them is less than one percent. However 
we cannot reproduce the fourth rank Hj moment of the data from our calculation. This 
result indicates that the oscillation of the Hj moment is very sensitive to the value of Fj 
moments. 

In order to see the effect of truncation for the normalized factorial moment Fj , those 



calculated with Eq.(p4|) (without truncation) and with truncated NBD are shown in 
Fig. 6b. The former is calculated with (n) = 8.80 and C^ = 1.190. The parameters of 
the truncated NBD are the same with those in Fig. 5. The Fj moment with truncation 
is much more suppressed than that without truncation at higher rank j. 

From Figs. 6a and 6b, we can see that the Fj moment calculated from the FFPE 
is smaller than that with Eq.(|34[) (without truncation). Therefore, introduction of the 
fractional derivative in time variable suppresses the value of Fj moment compared with 
that with Eq.(34) as the rank j increases, and gives rise to similar effect as in truncation 



of multiplicity distribution. 

The Hj moments and normalized factorial moments in ep collisions in the interval 
185 < W < 220 GeV ||l^ are also calculated with sets of parameters, a =0.25 and 
7t" =16.0 (A =19.94), or a =0.75 and 7t" =5.70 (A =20.41). The resullts become 
almost the same with those shown in Figs. 5 and 6. 

5 Concluding Remarks 

The fractional Fokker-Planck equation (FFPE) corresponding to the birth and death 
process with immigration is taken as a model for particle production processes with a 
memory effect. It is solved according to the procedure proposed by Barkai and Silbey [p| . 
From the solution of the FFPE, we obtain the generating function for the multiplicity dis- 
tribution, where parameter a connected with the fractional time derivative is contained. 
If Q is put to 1, the distribution becomes the NBD. 

The normalized factorial moment Fj calculated with Eq.(O) becomes much smaller 
than that with Eq.(^) obtained from the GF for the NBD as the rank j increases, where 
Eq.(pT|) coincides with Eq.(^) if a = 1. This fact means that the high multiplicity 
component in the multiplicity distribution is surppressed if a < 1, compared with the 
NBD (a = 1) with the same (n) and C^. 

When a is less than 1, the oscillation of Hj moment appears, and as a decreases 
from 1 to 0, the oscillation becomes much stronger. This is caused by the fact that the 
fractional derivative (0 < a < 1) is introduced into the time derivative in Eq.(y). It can 
be said that introducing the fractional derivative gives similar effect on the normalized 
factorial moment and the Hj moment as in the truncation of multiplicity distributions. 
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Figure captions 

Fig. 1 The Mittag-Leffler function calculated from Eq.(|2^) with a =0.25, 0.50, 0.75 and 
1.00. 

Fig. 2 Calculated results of Hj moments as a function of rank j] a) with a=0.25, and 
7f° =1.5 (A =48.44), 2.5 (16.93) and 3.5 (11.61). b) with q=0.50, and 7^° =1.0 
(A =66.82), 1.5 (15.41) and 2.0 (10.24). c) with q=0.75, and 7^" =0.5 (A =20.41), 
1.0 (10.20) and 1.5 (7.65). In each calculation, (n) = 29.2 and C^ = 1.274, observed 
values of charged particles in pp collisions at y/s =546 GeV/c ilj] are used 



Fig. 3 Calculated result with a = 0.5 and 7t° = 2.23 (A = 9.16) is compared with the 
data of charged particles in pp collisions at y/s =546 GeV/c [^], where (n) = 29.2 
and C^ = 1.274 are used. Parameter a is fixed at 0.5, and 7^" is adjusted so that 
the first relative minimum of the calculated Hj moment is located at j > 5 (j = 7) . 

Fig. 4 Calculated Hj moments with a=0.50 and 7t° =90.0 (A = 25.36), as a function 



of rank j are compared with those in e^e collisions at y/s = 91 GeV/c |15]. 
(n) = 20.70 and C^ = 11.091 are used in our calculation. 

Fig. 5 Calculated Hj moments as a function of rank j are compared with the data in 
ep collisions in the interval 185 < W < 220 GeV p!^ . White circles denote the 
calculated result with a=0.50 and "ft"' =11.2, where (n) = 8.80 and C^ = 1.190 
are used. White triangles show the Hj moment calculated with the NBD which is 
truncated at n = 21, the highest observed charged multiplicity, where (rib) = 8.825 
and A = 13.21 are used. 

Fig. 6 Calculated normalized factorial moments Fj as a function of rank j are compared 
with the data in ep collisions in the interval 185 <W< 220 GeV |l^. a) White 
circles denote the calculated result with a=0.50 and 7t" =11.2 (A = 20.00), and 
black circles are the data. 

b) White circles denote the result calculated with the NBD which is truncated at 
n = 21, where {n\^) = 8.825 and A = 13.21 are used. Crosses show calculated Hj 



moment from the GF, Eq.(15), where A = 13.10 is used. 
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